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A novel discrete growth model in 2+1 dimensions is presented in three equivalent formulations: 
i) directed motion of zigzags on a cylinder, ii) interacting interlaced TASEP layers, and iii) growing 
heap over 2D substrate with a restricted minimal local height gradient. We demonstrate that the 
coarse-grained behavior of this model is described by the two-dimensional Kardar-Parisi-Zhang 
equation. The coefhcients of different terms in this hydrodynamic equation can be derived from the 
steady state flow-density curve, the so called 'fundamental' diagram. A conjecture concerning the 
analytical form of this flow-density curve is presented and is verified numerically. 
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It is well established in the last two decades [H-Q that 
the one-dimensional Kardar-Parisi-Zhang (KPZ) equa- 
tion adequately describes the long-range dynamics of 
a collective motion of hopping particles on a line known 
as "Asymmetric Simple Exclusion Process" (ASEP) @- 
[^. Apart from many fruitful theoretical advantages, this 
ASEP-to-KPZ mapping enables a fast and simple way 
of modelling the KPZ dynamics. The latter is of wide 
interest since it appears in various contexts, (provided 
the symmetry of a nonequilibrium statistical system un- 
der discussion allows for the effective l-|-l-dimensional 
description), including, to name but a few, the models of 
crystal growth llOll . Molecular Beam Epitaxy . Burg- 
ers' turbulence I12l.[l3l. polynuclear growth p4l - [l8j . bal- 
listic deposition p^| - |22| etc. 

It is therefore an appealing idea to seek for a similar 
simple discrete multi-particle system whose long-range 
dynamics would be governed by a two-dimensional KPZ- 
type equation. Lately, several models of the desired na- 
ture, i.e. ones which combine the discreteness with the 
long-range KPZ-type dynamics, were suggested p3l - l28| . 
In this Letter we propose another model belonging to this 
class, which, in our opinion, combines the advantage of 
physical transparency with the flexibility of tuning the 
internal parameters of the model to catch the different 
desired regimes both in 1-1-1 and in 2+1 dimensions. This 
model, which we call "Zigzag Model" , has a simple geo- 
metrical formulation. 

Take an infinite cylinder covered by tilted square grid 
as shown in FiglT]a,,b and consider a directed closed path 
("zigzag") around a cylinder. Any such path consists 
of a constant number of rises and descents, constituting 
"kinks" (here and below we conventionally define rises 
and descents with respect to a rightward step). The den- 
sity of descents, p, is defined by the tilt angle a - see 
FiglT^,b where this density equals p =1/2 (for tana = 1) 
and p = 1/4 (for tana = 1/3), respectively. Consider 



evolution of a system of such nonintersecting zigzags. At 
each infinitesimal time step dt any elementary kink ori- 
ented downwards, can turn upwards with probability p dt 
under the condition that such a move is not blocked by 
the upper nearest neighboring zigzag (i.e. all zigzags stay 
nonintersecting at all times). By an appropriate rescal- 
ing of time, t, we set p = I. The examples of elementary 
jumps which are allowed and those which are blocked are 
shown in Fig|TJ 
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Figure 1: The zigzags on a cylinder for different tilt angles 
Of = 7r/4 (a) and a = arctan(l/3) (b), the examples of allowed 
an forbidden moves are shown. 

For better understanding of the dynamics of the model, 
notice that the evolution of a separated zigzag, (if, for 
the time being, we neglect its interaction with the other 
ones), can be interpreted as a hopping dynamics in a 
standard one-dimensional Totally Asymmetric Simple 
Exclusion Process (TASEP) ^ - see FigHfor the corre- 
sponding mapping, which is both conventional and self- 
explanatory. Essentially, a descent (rise) in a zigzag 
is identified with a particle (hole) in the corresponding 
TASEP. Therefore, the set of zigzags can be viewed as a 
system of interacting TASEP layers. 

The connection between zigzags and TASEP layers is 
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shown in Figl^Ji. Here, there are two sorts of constraints 
on the movement of the particles in layer B2 ■ Indeed, for 
a jump to be possible at some point of the A2 zigzag, two 
conditions should be simultaneously fulfilled: i) it should 
be a downward kink, and ii) the movement should not 
be blocked by the upper zigzag. The translation of the 
first condition into the TASEP language is convenient: 
there should be a particle at a given position in B2 and 
a void immediately to the right of it. In turn, the second 
rule (i.e., the one describing the interaction of the Bi and 
B2 layers) can be formulated as follows: it is possible to 
label the particles in the layers Bi and B2 in a way such 
that the fcth particle in the layer B2 can never surpass 
the particle with the same label (i.e., the k-th one) in the 
upper layer Bi, giving rise to an "interlaced TASEP" 
picture shown in Figl^jD [29j . 



Note, that Ay(x) > 1 at all positions. Therefore, for each 
particle on the upper line there is a "partner" particle on 
the lower line that cannot surpass it, and this "partner- 
ship" is preserved by the elementary moves. Indeed, if 
there is a particle on the upper line Ai in position x — 1/2 
(particle positions are shifted by half-step as compared 
to the positions of the kinks) and the zigzag-to-zigzag 
distance at x is A.y{x) then each particle on A2 jumping 
from .T — l/2toa; + l/2 will decrease Ay by one, and, 
therefore, up to (Ay — 1) particles can pass this point 
without being affected, while the particle number Ay, 
which is the desired "partner" particle, will get stacked 
at X — 1/2. Moreover, assume now that the upper parti- 
cle hops from a; — l/2to2; + l/2. For this move to take 
place, the position x + 1/2 should be a void before it, 
thus (compare to EqlT]) 
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Figure 2: The interconnections between zigzag, interacting 
TASEP and tag diffusion representations of the model, a) 
Interaction of two zigzags and the two corresponding TASEP 
layers, examples of the same allowed and forbidden moves are 
shown in both representations; b) interacting TASEP layers 
with examples of allowed moves the bold dashed lines connect 
the particles with the same number (see explanation in the 
text); c) the system of tagged particles corresponding to the 
TASEP layers depicted in figure (b), with the same allowed 
moves highlighted; d) the result of the shift hij — >■ hi,j + j. 



To prove these jumping rules, consider the distance 
Ay between the zigzags Ai and A2 (see Figl5^). At two 
adjacent spatial positions it satisfies 

' Ay(a; — 1) + 1 there is a particle at x — 1/2 
on B2 and a void on Bi 

Ay{x) = < Ay(x' — 1) — 1 there is a void at a; — 1/2 

on B2 and a particle on Bi 

, Ay(a; — 1) otherwise 

(1) 



Ay(a; + 1) = 



Ay(a;) + 1 there is a particle at a; -f 1/2 
on B2 

Ay{x) otherwise 

(2) 

and Ay{x + 1) does not change as a result of the move. 
In both cases of ([2]) the "partner" particle after the move 
is still the same: either the distance and the order of 
particles do not change, or the distance increases but the 
(Ay -f- l)th particle as viewed from the new position is 
the same as Ayth from the old position. 

In case of a finite cylinder with periodic boundary con- 
ditions the "partner" particle can, in principle, lag behind 
by a whole lap. In this case, one should ensure that the 
two particles do not interact even if formally they are at 
the same place, so one should keep track of the "real" dis- 
tances between the particles (i.e., those which correspond 
to the distances between the kinks on the original cylin- 
der), not the "apparent" modulo N distances (compare 
to P). 

There is now evidence of some spatial symmetry in the 
model: for each given particle there are exactly two other 
particles, one to the right and one (with the same num- 
ber) "on top" of it - see the dashed lines connecting the 
particles of the same number in Figl^b, which can block 
its movement via excluded-volume interaction. To bet- 
ter exploit this symmetry, it is convenient to reformulate 
the model, following the logic of the so-called "tagged 
particle diffusion" introduced in @ in order to show that 
on a coarse-grained level the ASEP dynamics is subject 
to the one-dimensional KPZ equation. 

Consider a set of m TASEP layers (zigzags) of length 
N with n particles within each layer. Note that instead 
of enumerating the cites of TASEP layers and marking 
which particular cites are filled with particles, one can 
store the very same information in a different way by 
enumerating the particles with two indices i S [l,?^], 
j € [Ijin], and ascribing to each particle a 'height' hij 
equal to its position on the corresponding layer (as mea- 
sured from some arbitrary chosen "1st cite"). It is clear 
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now that locally the values in the hij matrix are increas- 
ing in i direction and non-decreasing in j direction. To 
make the model completely symmetric make the transfor- 
mation (compare [3l|) hij hij+j which (see FigI5}:,d) 
ensures that hij is now increasing function in both di- 
rections, i and j. Under this transformation the dynamic 
rules become totally symmetric: 



hj 1 + 1 with Prob. dt, if 



hij otherwise 



(3) 

One can interpret these rules as the dynamics of a heap 
growing over a two-dimensional substrate, with an addi- 
tional constraint of heap gradient being not less than 1 
in each of the transverse directions. 

To obtain insights about the large scale dynamics it 
is useful to consider a coarse-grained hydrodynamic de- 
scription of the model [1]. To describe the model at a 
coarse-grained level we will henceforth use (x, y) as the 
spatial coordinates instead of in the lattice model. 
The coarse-grained dynamics can be described by a local 
'smooth' velocity field u{x^y,t) which is, in fact, noth- 
ing but the average (over space-time element dxdydt) 



rate of successful hops hi 



hij + 1, i.e. the average 



probability to find simultaneously {h. 



i+ij - hij) > 1 and 
> 1. This value obviously depends on 
the average slope (or equivalently on the local densities 
in both directions in the Zigzag model) of the surface 
in both directions. Similar to the one-dimensional case 
[3, 3 we assume that the velocity field u{x, y, t) depends 
on the space-time coordinates only through the local den- 
sities, i.e., 



Assuming now | | 1 and 1 one can expand 

both px^y and u{px,Py) as power scries in the derivatives 
of d, obtaining up to the second order: 



- / _ dd s2 f dd 

Px{x,y,t) ^ p^n^ Px— + {Px) ( — 



/ ,^ - I , - dd ,_ ^2 / dd\^ 



which when substituted into ([S]) gives 



(7) 



dd 



dd 



U{Px,Py) = U{Px,Py) ~ P^Ux{p.x,Py)—~ PyUy{p^,Py) — 



+ (yPlux{Pa:,Py) + ]^'ptuxx{Px,Py)^{^ 

+ {vljUy{{'Px,'Py) + ]^T'lUyy{-px,-py)^ ^^^^ 

-2—2 /— — sdd dd 
+ PxPyU^viPx.Py)-Q^-Q^ + - 

(8) 

where for brevity we introduced the notation Ux, Uy, etc. 
for the partial derivatives 



du 



This allows us to write down a time-dependent dif- 
ferential equation for d which belongs to the two- 
dimensional KPZ class: 



u{x, y, t) = u{px{x, y, <), py{x, y, t)) 



(4) 



where Px^{x,y,t) = ^^'^g;^'*^ and Py^{x,y,t) = ^^^g^^'*^ 
Note that the value px introduced here has a meaning of 
particle density in the corresponding TASEP layer. 

Suppose that the slopes px,y are weakly fiuctuating 
around some average values 'Px,yi which are determined 
by the boundary conditions. Introduce the displacement 
d{x, y, t) of a particle located at {x, y) from its average 
position at t = 0: 



d{x, y, t) = h{x, y,t) - ^ - ^ 

Px Py 



(5) 



Then, in an analogy with p33 | the local values of px,y and 
d are connected via 



1 



1 dd{x,y,t) 



1 



1 dd{x,y,t) 



Px{x,y,t) 



dx ' py{x,y,t) Px 



dy 



(6) 



If = DAd + r^ + u{px,Py)^DAd + ?] + u{px,Py)~ 

Edd \ - \ - dd dd 

a=x,y a=x,y p=x.y 

(9) 

where D is a diffusion coefficient, r]{x,y^t) is the noise 
term 



{vix,y,t)) =0, {r]{x,y,t)ri{x',y',t')) = 5x,x'5y.y'5t,t' 



and Cq,, bap for {a, j3} = x,y are: 



(10) 



Plua{Px,Py), 

PluaiPx^Py) + ^pluaaiPx,Py)i 
PlpllJ'xyiPx^Py) 



(11) 
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Figure 3: Comparison of the theoretical and numerical results for u{px,py). a.)u{px,Py) dependence as measured by direct 
numerical simulations (see main text for the details about the simulations); b) the discrepancy between the numerical result 
for u and the one conjectured in Eq. 1151 c) the u{p) behavior in the p = px = Py plane, the points correspond to numerical 
results, while the dashed line is the naive mean-field guess u = (1 — p)^ (compare Eq. I13[) and bold line is the conjectured 
dependence u = (1 — p)(l — p/2) (compare Eq. I15|) : . 



The behavior of the system is controlled by the usual 
flow-density dependence u{px, Py), which is central in all 
models of traffic (see, for example, [s^) and is usually 
referred to as the 'fundamental' diagram. Indeed, the 
flow of the particles in the interlaced TASEP formulation 
of the zigzag model is equal to I{px,Py) = Pxu{px, Py)- 
Recall that in the ID-case the function u{p) = I{p)/p is 
just u = 1 — p and it can be obtained by the mean-field 
arguments. In absence of interaction between layers one 
would have expected the same dependence in our model: 

u{Px,Pv) = - Px (12) 

In the presence of interaction the corresponding 2D 
mean-field result would be 

u{Px,Py) = {'^ - Px){l - Py) (13) 

where the symmetry of the model is taken into account: 
a particle can hop if there is a void both in front (with 
probability 1 — p^;) and on top (with probability 1 — py) 
of it and the horizontal and vertical jumps are supposed 
to be independent. However, contrary to the ID case, 
this result is not exact. Indeed, the connectivity of the 
surface dictates that the local increments (/ii+i,j — h.i j-^ 
and (/iij-|_i — hij) are positively correlated and cannot 
be considered as independent. 

The exact analytical evaluation of the velocity 
u(px,Py) up to now is beyond our reach. We have made 
numerical simulations of u{px, Py), the results being pre- 
sented in Fig [31 The simulations were performed for the 
systems of size iV = 32 in both directions with periodic 
boundary conditions ensuring that the average densities 
take the values Px = 7>y = 0, 1 [l^l- The results 

were averaged over 9 x 10^ Monte-Carlo steps. To make 
the comparison with the mean field more visually com- 
pelling, we plot in Fig|3t separately the numerical data 
for Px ^ Py = p together with the mean-field results (fT2)) 
and p^ . As expected, the first of them overestimates 



the flow, while the second one underestimates it. In fact, 
the numerical data fits perfectly the form 

^.(p) = (l-p)(l-p/2) (14) 

and the consideration of limiting cases p — > and p — >■ 1 
suggests that this result may be exact. In particular by 
developing the perturbation theory at high densities, we 
are able prove that u{p) -> (1 — p)/2 as p — > 1, the details 
of these computations will be published separately [35i ] . 

The simplest generalization of Eq. p^ onto the case of 
Px 7^ Py which respects the boundary conditions it(0, 0) = 
1, m(1, £) = u(^, 1) = for any ^ S [0, 1] reads 

(15) 

In the FigOjD we have plotted the discrepancy 
Au{px,Py) = Ua{px,Py) - u{px,Py) between numerical, 
UniPxiPy) and conjectured, u[px,Py), (see ([T5|) ) func- 
tions. One sees that u{px, Py) is in very good agreement 
with the results of numerical simulations. 

The conjectured function u{px,Py) allows to evaluate 
the coefficients in (|lip for any px and py and to cal- 
culate the eigenvalues of the coefficients matrix hap in 
front of the nonlinear term in Eq.(jni). We show thus 
(see FigH]) that the domain [0 < p^^ < 1, < pj, < 1] 
can be separated into two regions: i) a region where one 
of the eigenvalues dominates (|Ai| ^ IA2I) signahng a 
quasi-one-dimensional behavior, and ii) a region where 
I All ^ IA2I corresponding to the truly two-dimensional 
KPZ dynamics. In the latter region both eigenvalues are 
negative, and thus the matrix is positively defined. 

Summing up, we have presented a novel model of a 
statistical driven system in 2-1-1 dimensions which turned 
out to be a direct generalization of the conventional one 
dimensional TASEP model. In the hydrodynamic limit 
we derived the differential equation for the particle dis- 
placement in this model, and we have made a conjecture 
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Figure 4: Density plot of the ratio A1/A2 of the eigenvalues 
of the matrix B. In the region det B < the eigenvalues have 
opposite signs. Outside the dashed area the system can be 
regarded as effectively one-dimensional since Ai ^ A2. 



concerning the flow-on-density dependence of the model 
(the so called 'fundamental' diagram). Clearly, the model 
suggested here deserves further investigations in various 
directions. In particular, a detailed investigation of the 
limiting cases Px = Py —> and p^, = pj, — >■ 1 can give 
valuable understanding of the flow-density dependence, 
the study of the displacement fluctuations in different 
regimes (compare to |3|) should be very fruitful, and the 
problem of proving conjecture or at least sounds 
quite challenging. Some progresses in these direction will 
be reported in a forthcoming longer paper (ssl). 
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